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Abstract 

We numerically compute features in the power-spectrum that originate from the decay of helds during 
inflation. Using a simple, phenomenological, multi-field setup, we increase the number of fields from a few 
to thousands. Whenever a field decays, its associated potential energy is transferred into radiation, causing 
a jump in the equation of state parameter and mode mixing at the perturbed level. We observe discrete 
steps in the power-spectrum if the number of fields is low, in agreement with analytic arguments in the 
literature. These features become increasingly smeared out once many fields decay within a given Hubble 
time. In this regime we confirm the validity of the analytic approach to staggered infiation, which is based 
on a coarse-graining procedure. Our numerical approach bridges the aforementioned analytic treatments, 
and can be used in more complicated scenarios. 
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INTRODUCTION 



flfl 



The construction of inflationary models in string theory is an active research field, see [iHi for 
reviews. Early proposals often involved few dynamical degrees of freedom, i.e. a single field dri ving 
inflation, while all others are meticulously stabilized as in the KKLMMT [tI brane inflation js- 15| 
setup. These models offer computational control and predictability, but they appear simplistic in 
the absence of any a priory reason or need for such frugality; for instance, if inflation is driven 
by a brane/anti-brane pair, why shouldn't mo re p airs be included? As a consequence, multi-field 



models have become increasingly popular 



24 



especially since observable non-Gaussianities, 



35l |. are possible [3g] (see also 



for which there is emerging experimental evidence 

If inflation is driven by more than one degree of freedom, the end of inflation can differ signifi- 
cantly: instead of a sudden end caused by the simultaneous decay of all fields, a stretched out decay 
38 h40 | , a nd in some cases unavoidable and the feasibility of (p)reheating 



phase is possible 



can change drastically 



41 



4J]. Staggered inflation 



3 



is not a new type of inflation, but a 



collective term coined for models that contain such decaying fields during inflation. If all fields 
decay in a few e-folds, as in [28], the preheating phase is merely extended, without strong obser- 
vational consequences, but if they decay throughout the last sixty e-folds, additional signatures in 
the correlation functions of fluctuations in the Cosmic Microwave Background (CMB) are possible 



38h4CII| (see also j45l447l| for additional signatures caused by particle production during inflation. 



as present in trapped inflation 



mIsqI]). 



We use the term "decay" loosely to indicate the partial or full transfer of an inflaton's energy 
(potential and/or kinetic) to an additional component of the energy momentum tensor, such as 
radiation. A concrete example is the annihilation of a brane/antibrane pair in the extension of 
the KKLMMT proposal in j33[; whenever a brane comes in close contact to an anti-brane, they 
annihilate and produce closed string modes that redshift like radiation. Other examples include 
inflation driven by tachyons 
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32( 1 . by multiple M5-branes [28h3C1I| or inflation on the landscape 
[40], among others. 

If many fields decay in any given Hubble time, an analytic formalism, based on coarse graining, 
was developed in [3S 



40l ]: the aim of this study was to retain some of the effects caused by the 
decaying fields, such as contributions to the power-spectrum generated by the additional decrease of 
the energy that drives inflation, due to the fields' decay as opposed to slow roll. As a consequence, 
contributions proportional to T/H, where T is the decay rate and H the Hubble parameter, appear 
alongside the common slow roll parameters in observables such as the scalar spectral index or the 
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tensor to scalar ratio; furthermore, the new contributions may even be the dominant ones 40l |. 

However, any effect due to the sudden decay of an individual field is not retained. To ameliorate 
this shortcoming, a single decay was discussed in detail in a concrete setup (an extension of the 



3j]. There it was shown that 



KKLMMT proposal) in [33] and in a related DBI-inflation setup 
a jump in the equation of state parameter results after the decay. This discontinuity in turn 
causes mode mixing at the perturbed level ^. In 33|] the matching conditions for perturbations 



were derived in the sudden decay approximation, neglecting perturbations in radiation. A ringing 



pattern on top of the nearly scale invariant power-spectrum was found for sub- horizon modes [33 1 
(see also [5^ for the same effect induced not by decays, but by a var ying speed of sound) and super- 
horizon modes showed small corrections. Extending the model of |33l] f urther, to allow for many 
decays, should yield a power-spectrum that converges to the one of [38|, [4^ even for sub-horizon 
modes, as a superposition washes out any signals of the individual decays. 

To this end, we develop the numerical tools needed for such a comparison jn this paper. We 
apply them to a much simpler setu p m otivated by inflation on the landscape 
a case study for the formalism in 



40|, that served as 



331 ] in the near 



and enables us to compare results. We plan to use the 
code developed for this study in the more complicated, yet more realistic model of 
future. 

Besides developing and testing the numerical code, we show how the analytic treatment of 
staggered inflation in 



3S 



401 ] is recovered in the large AA-limit, and argue that this formalism 



serves as an excellent approximation if more than ~ 40 fields decay in any given Hubble time. 



The outline of this paper is as follows: in Sec. Ill Al we specify the mode 



analytic results of 



3g 



401 ] in Sec. IIIBI The matching conditions derived in 



before reviewing the 



33(, which are needed 



for the numerical treatment, are provided in Sec. lIH Al We then show in Sec. lIVI how the relatively 



simple results of 



401 ] are recovered numerically in the large A/"- limit, while we recover signatures of 



individual decays if only a few fields are present. Details about the numerical code can be found 
in the Appendix. 

If not stated otherwise we set the reduced Planck mass equal to one, m~j^ = SttG = 1 and use 
"~" to denote equality to leading order in small parameters. 



^ Mode-mixing means that the second independent solution of the Sasaki-Mukhanov variable, which has zero am- 
plitude if Bunch-Davies vacuum initial conditions are imposed, acquires a non-zero amplitude. 
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II. STAGGERED MULTI-FIELD INFLATION 



A. The Background Model 

Consider M uncoupled, scalar fields (f)i, I = l...J\f with canonical kinetic terms and linear 
potentials 

Vi = Vq' - ci^i . (1) 



Such a potential is motivated 



that reside on flat stretches 



by expanding a general potential on the string-landscape for fields 



33( 1 : fields located on steep slopes evolve faster and become dynamically 
irrelevant for the evolution of the universe. Hence, we are interested in a narrow distribution of 
the cj and fields with comparable potential energies. To keep the model as simple as possible, we 
consider cj = cj = c and Vq = Vq for all /, J (see [ss] for more general setups). 

Such an expansion is expected to be valid for small ranges of field values only, so we use ([1]) 
only up until some mclxiinal field, value (pj — 4^end 

^ 1; if a field encounters (j)end, we assume that 
its potential energy is converted into an additional component of the energy momentum tensor, 
i.e. radiation pr. Thereafter, the field may be stabilized or continue to roll freely. In either case the 
field quickly becomes irrelevant for the evolution of the universe. We choose to set the potential 
to zero for 0/ > (pend ^-ncl let fields roll freely for larger field values, always keeping (pj and (pj 
continuous through time^. Since we would like the previous stretches to be flat, we demand 

(2) 

Depending on initial conditions, some fields encounter (pend earlier than others, causing "decaying" 
or "dropping" fields to enter the freely rolling phase in a staggered fashion. As we have no knowledge 
about the concrete initial conditions, we distribute the initial field values randomly over the interval 
{0,(pend)j from which they evolve according to their Klein-Gordon equations 

dVr 

cPj + 3H(Pj = -—f, (3) 
o(pi 

with V from ([1]) for (pj < (pend and ^ = otherwise. Radiation is produced whenever a field 
encounters (pend-, taking over the remaining potential energy and redshifting as 

Pr = -AHpr (4) 

thereafter. 



^ Since the kinetic energy is negligible compared to the potential one, we do not expect any significant changes if 
paxt or all of the kinetic energy were also infused into radiation. 
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B. Analytic Results in the Large M Limit 



If many fields drop out during any given Hubble time, one can describe the system analytically 
jssl . I40I by promoting the number of fields to a continuous function that decreases in time according 



to the decay rate 



r^-^>o. (5) 



This rate is not a new free function but set by the initial values and subsequent evolution of the 
fields. Once T ^ H, inflation ends because all remaining fields decay in a single Hubble time. 
Hence, inflation driven by many fields can take place as long as T <^ H , that is as long as [38] 

eM = ^<^l, (6) 

where we included a 2 in the definition of ej\f to simplify expressions below. 

If inflation is of the slow roll type, that is if the slow roll parameters £/ = ( dVT /dtpr)'^ /(2W'^), 



5l| for details on the 



e = ^^^^ VI = (d'^Vi / dip'j) /W and rj = J^i'Hi small (see 

Hubble slow roll expansion for multi- field inflation), one can show that the combined energy 
of the fields and the one of radiation obey 

Anf ^ -2H{e^ + e)pi^t, (7) 
Pr ~ -2H 0(1 + Wr)pr - £j\fPin^ ^ '^H{ej^ - e)pi,,t , (8) 

to first order in small parameters. Here Wr = Pr/ Pr = 1/3 and 

e=hl + wr) ^2 ^EN- (9) 

Z Pr + Pinf 

In the limit of many decays during any given Hubble time, radiation approaches a scaling solution 
where the energy loss due to redshifting is compensated by the infusion of energy from the inflaton 
sector, Pr ~ e]^2pij(?, + 3^1),.). Further, the Hubble slow evolution parameter becomes 

-^o,e + e. (10) 

The analysis of adiabatic perturbations is straightforward, albeit tedious, and one can show that 
in certain cases (for example in the model we investigate in this paper) the effects of isocurvature 



38| 



perturbations are negligible [38]. T 



le scalar power-spectrum of the curvature perturbation on 



uniform density surfaces Ck becomes 



87r2m2;(e72 + e) ' ^^^^ 



where we reinstated^the reduced Planck mass and 7 is a parameter of order one set by the back- 
ground evolution [38]; since we are primarily interested in situations where the staggered inflation 
effects dominates over slow roll effects, e,r/ <C e, we do not need to compute 7 and the scalar 
spectral index in 



381, 



40| simplifies to 



where 



TH 

Similarly, the tensor power-spectrum can be computed to 



~ (5-3)e, (13) 



(14) 



VT^I^AhW (15) 



■(2^)3 



with the tensor spectral index 



and the tensor to scalar ratio 



For the linear potential in ([T|), we get 



SO that 



(16) 



~ -2e. (18) 



r = ^ (19) 
~ 16e. (20) 



r^ ., ".n (21) 



'^^^""-la^i^Zv^H- ^^^^ 

and 

(5~-l. (23) 

The number of e-folds ~ 60 becomes 

iV = r'i^dt:.^^, (24) 



which allows us to express e in terms of N, 

The resulting scalar spectral index (fTSll . tensor spectral index (fTSjl and tensor to scalar ratio (f20l) 
are 

(26) 

nrc.-^, (27) 
r.l. (28) 

These expressions are valid if the slow roll contributions are subdominant and many fields drop out 
in any given Hubble time so that the coarse grained analytic treatment is justified. All signatures 
due to the sudden changes in the equation of state parameter whenever a field drops out are not 
retained. 

C. Model Parameters 

In the next section we aim to test these analytic predictions numerically, without performing 
a coarse graining. To this end, we need to specify the model parameters. To avoid confusion, 
we restore the reduced Planck mass rupi = (SttG)^^^"^ in this section. First, we want to compute 
observables around sixty e-folds before the end of inflation when modes that are observable in the 
CMB crossed the Hubble horizon, 

iV = 60 . (29) 

We will vary the number of fields J\f, but are particularly interested in the limit J\f ^ 1. Further- 
more, we would like to impose e, c/(/)end/M) ^; to be concrete we set 



-pi \/^7VU '-Yij/ ^""pi 

SO that 



17 ' (30) 



The inflationary scale follows from (|lip with the COBE normalization V( ~ V-r. ~ 2.4 x 10^^, so 
that 

^0 = 3^2.4 xlO-^m^,, (32) 



where we used e ^ e. The critical field value follows from (j24p as 

1 



'fend 



2VM 



m. 



pi ' 



(33) 



One may check that c/(/)end/Vb = e/{2\/JJ) <^ e, as desired. 



III. PERTURBATIONS 

We now go beyond the analytic approximations in jssl. |40|| . such as the large AA-limit and the 
slow roll approximation, and investigate the evolution of scalar perturbations numerically. 
The line element to linear order in scalar perturbations and without fixing a gauge is 



ds"^ = -(1 + 2A)dt^ + 2aB^idx'dt + a'^[{l - 2i))5ij + 2E^ij]dx'dx^ 



(34) 



Metric degrees of freedom couple to the perturbations of the scalar fields 6(j)j in their equations of 



motion 



53, 



bj + 3H6(j)j + —d(j)i + V<i>j<t>j5<l)j = -2V^^A + 0/ 



k' 

A + Sip + ^l a^E - aB 



(35) 



with / = l,...,J\f. Two metric degrees of freedom can be eliminated by the choice of gauge. 



Utilizing this choice, we use two gauge-invariant perturbations, i.e. the two Bardeen potentials 



55| 



^ = A + [aB - a^Ej ' , 
^^^ = ^P-H (aB - a'^s) . 



(36) 
(37) 

Since anisotropic stress is absent in our setup, we have $ = ^. It is also useful to introduce gauge 
invariant field perturbations, for instance the Sasaki-Mukhanov variables, 



Qi 



(38) 



where t 
Qj are 



re spatially flat gauge is defined by the condition ip = 0. The equations of motion for the 



53 



63| 



= Qi + 3HQi + —Qi + YIv 



<t>i4>j 



' TT'^I'^J 1 1 Qj 



mlia^ V H 



(39) 



in the absence of perturbations in additional degrees of freedom. Even though we keep track of 
radiation at the background level, we ignore perturbations in pr. 
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The Sasaki-Mukhanov variables have the advantage that their evolution is decoupled from the 
metric degrees of freedom. We assume that the Qi are in the Bunch Davies vacuum at the onset 
of inflation, for —kr — >■ oo, that is 

QKtt) = ^ = -^^e, . (40) 

a{Ti) a{Ti)V2k 

where the e/ are independent unit Gaussian random fields with 

<e, >=0, <ei{k)ej{k')>=6ij6^{k-k'). (41) 

Since our ultimate goal is to numerically solve the equations of motion (|39|) . we want to use analytic 
solutions up until the wavelength of a mode approaches the horizon size. However, whenever a 
field drops out, the analytic solutions do not conform to the vacuum solution any more, but carry 
an admixture of the second independent solution of ([5^ . 

Qiitt) = i^e^-(^+V^)/^ («<)(-fcr) + /3<)(-A;r)) e, , (42) 
^ ^aie-'^^ + /3ie'^A e, , (43) 



aV2k 

where we expand the Hankel functions for large arguments and used /i = 3/2 + 0{e,£,r]) ~ 3/2 
during inflation. The Bogoliubov coefficients a and /? need to be computed according to the 
matching conditions of the next section whenever a field drops. 



A. Matching Conditions 



le equation of state parameter makes 



331] where the matching conditions of 



When fields encounter the sharp drop in the potential, t : 
a jump due to the creation of radiation. We follow closely 
perturbations were derived in a related setup ^. Based on the continuity of the induced metric and 
the extrinsic curvatures on the hyper-surface at which the equation of state parameter jumps, one 



can derive the matching conditions for the Sasaki-Mukhanov variables [56H58I]. The transitions, 
which we model as instantaneous events^, occur at well defined values of the infiatons. As a result, 
the Bardeen potential and the comoving curvature perturbation are both continuous^ 59146 1[| 



(44) 



^ The simpler matching conditions in [3j| lead to jumps in the extrinsic curvature and/or induced metric on the 
transition hypersurface in contradiction to [stI-IHsI]. but might still provide a viable approximation for super-horizon 
modes. 

In more realistic scenarios, our treatment remains reliable for k <C {At)~^ where At is the time-scale of the 
transition. 

^ We ignore terms of order C'(fc^), since we are interested in the effects of inflaton decays during the last sixty e-folds 
onto super-horizon modes. 
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where 



7^ = $ + 



(45) 



3(i + w;) yn 

The subscripts — and + refer to the value of the quantity in brackets before and after the field 
drops out, respectively. At the background level, the scale-factor, Hubble parameter Ti = a' /a and 
velocities of the fields are continuous 



[a]± = o , [n]± = o , Wj]± = 0. 

The comoving curvature perturbation, TZ, is related to the Sasaki-Mukhanov variables via 

1 



(46) 



n 



(47) 



3(1 + w)ml^'H ^ 

If no additional contributions to the energy momentum tensor are present at the perturbed level 



we still need to express $ in terms of t 
need the perturbed Einstein equations 



le Qj to derive their matching conditions; to this end, we 



55| 



2ml ° 



1 



with the gauge invariant perturbations of the energy momentum tensor 



(gi) 

i ' 



and the gauge invariant field perturbation 



i^igi) 



bi + (t)'i{B - E') . 



After some algebra 



33] we arrive at 



(48) 
(49) 

(50) 
(51) 

(52) 
(53) 



In order to guarantee the continuity of TZ and "I?, we impose matching conditions for all elements 
of the sums in ()47p and (|53p individually, 

Qi 



l + w 



Q'i4>'i + Qi {a\4„ + w)(t>'j7?j 



0, 



0. 



(54) 
(55) 
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For ease of notation we define 

(I54p leads to an enhancement of Qj on all scales each time an inflaton encounters a critical value. 
The magnitude of the jump is determined by the ratio 

B^l±^. (57) 

If e ^ e, as in Sec JII CI the equation of state parameter becomes to leading order in e 

l + wc^le = 2 ^ , (58) 

3 Pr+ Pint 



so that 



(59) 

Pr ^ 



where AA^ ~ 0{N/N) is the number of e-folds between successive drop of fields. Thus, if 



many 



33], it 



fields drop out in any given Hubble time, M ^ N, the jumps become exceedingly small. In 
was found in a related model that the deviation of B from one also determines the amplitude of 
additional features in the power-spectrum. Based on this, we expect a requirement of N/N <C 1, 
that is at least M ~ 0(10^) fields, to remain observationally viable. 

Let's turn our attention to the second matching condition, (j54p . If the system is in a slow roll 
regime, the speed Q'j is continuous to zeroth order in small parameters for all but the one field 
whose potential energy has been infused into radiation. To see this, we insert the slow roll equation 
of motion, 2>}i(f)' = —a^V^j, into (j56p and work to zeroth order in slow roll to arrive at 

3 

Aii^-il + w). (60) 
Thus, using (I58p and working to zeroth order in e and the slow roll parameters, one may use 

[Q/]±^0, (61) 



instead of (j44p . This matching condition is a good approximation for the fields responsible for 
reheating the universe after inflation, which do not encounter (pend until the end of inflation^. 
These fields determine the scalar power-spectrum and are thus the ones of primary interest to us. 



^ We use the full matching conditions in our numerical code. 
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1. Matching the Bogoliubov Coefficients 



As long as — A;r S> 1, we rely on the analytic approximation for the Qj in (|43|) . Initially, aj = 1 
and f3j = 0, but /3/ 7^ as soon as the first field encounters the step in the potential. At each of 
those instances we have to match the analytic solution in (j43p according to (j54p and ()44p . yielding 



aT = a. - IH ^ + t— ^(1 - ^7) - (I- Ar) 



with A from (1561) and we used —tc^T-L . One can check that 



I +|2 |„+|2 1+'U;+/| _|2 |o-|2\ 



as in [33]. Note that all modes need to be renormalized by l/V B after the matching to guarantee 

\ajf -\l3+f = 1. 

There is, however, another subtlety related to the sudden decay approximation At^ — t- 0: in 
this limit all modes within the horizon are affected by the transition, leading to unsuppressed 
oscillatory corrections on top of the power-spectrum for all —kr ^ 1 [33]. However, in a realistic 
setting the infusion of potential energy from an inflationary field into radiation needs time. We do 
not expect modes with wavelength much smaller than this time-scale {k ^ (Ard)^^) to be affected. 
This in turn should lead to damping of the oscillations on top of the power-spectrum. Since we 
are working from a phenomenological viewpoint, we treat Ar^^ as a free parameter^. 

To incorporate this effect into the matching conditions (162]) and ([63]) we define a simple smooth 
window function 



F{k) = ^{1 + tanh 



c ( log ( A ) _ 1 



(65) 



The value of c determines the sharpness of the window function in k-space and kt indicates the 
limiting wave-number above which modes are unaffected by the decay. 



The source of \cef\^ — \l^t\^ 7^ 1 is our neglection of perturbations in radiation. 

In the extension of the KKLMMT model in [s^l , this time-scale is set by the brane decay time. 
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This window function can be used to modify the matching conditions to 



at = y^Cj[VF+\/l^aj \ , (66) 



/3/ = ' ^^^^ 

where we denote with af and f3f the right hand side of ()62p and ()63p respectively and we define 

Ci = ^ — -— . (68) 

l + ^F{l-F)/B {af + ap) 

The Ci and the factors of B guarantee the proper normahzation of the Bogohubov coefficients, 
|a^P — = 1- Note that the C/ only differ from one in a region around kt, whose width is set 

by the sharpness of the window function, c. We won't make use of the smooth window-function in 
our code, but use a sharp theta-function for simplicity. 



B. The Power-spectrum 



The power-spectrum of the curvature perturbation^ is defined as (see [6^ for a review 

W) 

and the scalar spectral index is given by 



53(k - ]i!)Vn = j—^ < 7^(k')*7^(k) > (69) 



ns-l = dln'PTz/d In k . (70) 



IV. NUMERICAL RESULTS AND COMPARISON TO ANALYTICS 

We solve the staggered inflation model of section III Al numerically with the model parameters 
of section III CI and a varying number of fields M = 1 . . . 2500. Perturbations are treated analyti- 
cally deep inside the horizon, according to (j43p . with or without adjusted Bogoliubov coefficients 
according to (j66p and ()67p whenever a field drops/decays. Not imposing the matching corresponds 
to a sharp window function at k^un^ = e^Hmi] once a mode approaches the horizon, we solve (f39|) 
numerically, together with the matching conditions in (I54p and (I44p . More details on the code can 
be found in appendix lAl Towards the end of inflation, we evaluate the power-spectrum in ()69p . 
which is dominated by perturbations in the fields that drove inflation the longest, and read off the 
spectral index. 

^ TZk coincides with C,k on large scales, the difference being terms of order . 
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How many fields need to decay in any given Hubble time for the large limit of Sec. IHBI 



401 1 to be a good approximation? To get a simple order of magnitude estimate, consider 
e = 2pr/{pr + Pinf) in & and assume that around QO/Af fields drop out in a Hubble time, e ~ 
C(10~^) needs to be smooth to at least one part in ten in order to describe properly the deviation 
of the power-spectrum from scale invariance, (j26|) . Lets focus on the relevant regime around 
sixty e-folds before the end of inflation when M is large. Since Apr = Pmf/-^ we get from 
Ae/e ~ 2/{MW'^) < 10"^ a lower bound of TV > 2000. Another way of getting an estimate is 
to consider how much pr redshifts until the next field drops out: since e ~ W-^{ao/a{t))'^ and 
a{t)/ao « ex.p{HAt) ^ 1 + HAt we get by using At ~ 60/7V a lower bound of TV > 2400. Thus 
we need around 

Ar>0{W^), (71) 

fields initially. 

We begin by simulating a small number of fields, J\f = 5. Fig. [1] shows the main quantities of 
interest from this simulation. 

In Fig. [Te] and Fig. [TJ we compare the effect of the window function onto the power-spectrum. 
In Fig. [Te] the Bogoliubov coefficients are fixed and in Fig. [lf| they are matched every time a field 
drops out according to ([66|) and ([67j) . Fixed Bogoliubov coefficients correspond to a sharp window 
function at k^um- Without this window function, all modes within the horizon are affected by 
a field decay (an unphysical effect), which explains the oscillations in Fig. [ifl If the number of 
fields is large, slow-roll is hardly disturbed whenever a field encounters (j)end and whether or not a 
window function is used becomes irrelevant. We fix the Bogoliubov coefficients for the rest of our 
simulations-utilizing a sharp window function. 

In Figl2] we plot e over the number of efolds for varying and in Fig. [3] we plot the Hubble 
parameter, H, for the same set of simulations. For small M, we observe a step-like reduction of H 
whenever a field encounters a drop. These steps are smoothed out, since we track the radiation into 
which the fields' potential energy is infused. This goes hand in hand with large changes in e. Fig. [2l 
The presence of these steps, and the accompanying change in the equation of state parameter, give 
rise to ringing patterns in the power-spectrum for large k. For low J\f, these ringing patterns 
are problematic due to their large amplitude. Away from the decay events, the power-spectrum 
becomes smooth again, as expected during slow roll. As the number of fields is increased, the steps 
begin to overlap, leading to a smooth reduction of the inflationary energy and less variations in 
e; as a consequence, features in the power-spectrum overlap, leading to a smooth spectrum in the 
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(a) 




(c) 
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(d) 
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(e) 



(f) 



FIG. 1: The first four panels are (a) tfie evolution of each (pj, (b) the Hubble parameter, 
H/rrirmpU (c) the equation of state \ + w and (d) e as a function of time throughout the 
simulation for M = 5. The final two panels show the power spectrum when (e) the Bogoliubov 
coefficients are fixed and (f) Bogoliubov are matched every time a field drops out. The simulation 
begins with ao = 1 and the power spectra are normalized to unity when k^, is horizon sized 60 

e-folds before the end of inflation. 



large J\f limit, see Fig. H] and Fig. [5l 

A word of caution is in order here: in our plots, we rescaled the amplitude of the power-spectrum 
to one at sixty efolds before the end of inflation. The actual amplitude increased in our numerical 
code (independent of k) as the number of flelds increases. This effect is due to the rescaling of the 
Bogoliubov coefficients whenever a field drops out and thus an artifact of ignoring perturbations 
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FIG. 2: The time evolution of e = 2pr/{pr + Pini) is plotted, using the number of e-folds as a time 

variable. The number of fields is varied, a) TV = 10, (b) Af = 25, (c) TV = 50, (d) TV = 100, 
(e)TV = 250, (f) TV = 500. For small TV, e changes drastically. For large TV, radiation approaches 
a scaling regime and the validity range of the analytic approach to staggered inflation in [s^, 4^ 

is entered. 



in radiation. The shape of the spectrum is entirely unaffected. 



is easy to understand, since the 



The presence of ringing patterns at large k as found in 
power-spectrum is modulated by the combination \a — /3p of the Bogoliubov coefficients. If we 
consider the case of a single decay (fSj = and aj = 1), work in zeroth order of small parameters 
{Aj — )• 0), and ignore the window function in ([65|) we arrive at a modulation of the I'th field's 
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FIG. 3: The time evolution of the Hubble parameter is plotted, using the number of e-folds as a 
time variable. The number of fields is varied a) M = 10, (b) M = 25, (c) M = 50, (d) M = 100, 
(e)AA = 250, (f) J\f = 500. For small J\f, step-like features are evident, which wash out once J\f is 
increased. For large M, the validity range of the analytic approach to staggered inflation in 
Qj 3] is approached as the evolution of the Hubble parameter becomes smooth. 



power-spectrum of 



I |2 



1 B-l f , ^ 1 

1 1 — cos(2x) sin(2x) 

B B \ X 

+ ^ ' ( 1 - cos(2x) ( 1 + — 



(72) 

2B V 

where we used ([62]) and ([63|) . rescaled the Bogoliubov coefficients with \l\fB to guarantee |af | — 
l/S/p = 1, and defined x = fe/^decay in terms of the comoving mavenumber k and the Hubble scale 
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(a) 



(b) 




(c) 



(d) 




(e) 



(f) 



FIG. 4: The power-spectrum evaluated numerically towards the end of mflation is plotted over 
the wave-number for (a) N = 10, (b) N = 25, (c) N = 50, (d) N = 100, {e)N = 250, (f) 
M = 500. Signatures of individual decays wash out with an increasing number of fields. 



at the time of the transition T^dccay Only modes within the horizon at the time of the transition 
pick up these modulations; if more fields subsequently decay, the patterns get more complicated 
and ultimately wash out. Furthermore, the oscillations are damped away once k approaches the 
inverse of the decay time. 

Another word of caution: since we used [7^]^ = in deriving the above, we ignored a term 
proportional to = (k'H)'^; if this term were kept, additional contributions would result for large 
X that scale as x^, dominating over the corrections we kept in the large k limit. Thus we provide 
a lower bound on features in the power-spectrum by using [TZ]^ = and the resulting matching 
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conditions. 

For7\A~ 10^, we start to recover the analytic results of jssl. kd] . as evident by Fig. [Hand [S] and 
a comparison of the scalar spectral index found numerically to ()26p . see table [H To compute the 
spectral index in tableHJ we average the slope of the power-spectrum in Fig. [5] and [5] between 60 and 
30 e-folds before the end of inflation. Naturally, such a smoothing is not a good approximation for 
the spectrum if pronounced patterns are present, as in Fig. [T]and Fig. [2] for low N . Thus the low 
M values for Ug should be taken with caution. Since the running time of our code scales linearly 
with A/", and we track several hundred modes per field^^, we were currently limited to ~ 10'^ fields 
(a run takes around a week). 



N 


Us - I 


5 


-0.0149 


10 


-0.0269 


25 


-0.0540 


50 


-0.0375 


100 


-0.0472 


250 


-0.0388 


500 


-0.0270 


1000 


-0.0165 


1500 


-0.00676 


2000 


-0.00598 


2500 


-0.00707 


oo 


-0.0167 



TABLE I: Comparison of the scalar spectral indices found numerically for varying number of 
fields to the asymptotic M = oo result in (f26]l . ris — 1 ~ — 1/iV = —1/60 —0.0167, based on the 

analytic treatment of staggered inflation in Sec. Ill Bl [s^, [4^. For low M features in the 
power-spectrum are so pronounced that a power-law is not a good approximation and the values 

for Ug should be taken with caution. The difference in the large N limit is caused by our 
neglection of perturbations in in our code, which are accounted for in the analytic approach. 

We see that the analytic approach to staggered inflation in [s^, becomes a good approxima- 
tion if around 40 or more fields drop out in a given Hubble time. In this case the evolution of H 



200 modes are tracked in the simulations in Fig. 2] and 50 modes are tracked in the simulations in Fig. El 
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(a) 



(b) 




(c) 



(d) 



FIG. 5: The power-spectrum evaluated numerically towards the end of inflation is plotted over 
the wave-number for (a) N = 1000, (b) N = 1500, (c) N = 2000, (d) N = 2500. These 
simulations involve considerably fewer modes than the models depicted in Fig. [3 



is sufficiently smooth when modes relevant to observations leave the horizon and features due to 
individual decays in the power-spectrum are washed out. 

To summarize, from our comparisons we can draw two conclusions; firstly, the analytic results 
of js^ . [40I can be recovered numerically in the large AA-limit. Furthermore, the distinct features 
for low N are in , our numerical code reliably tracks the evolution of 

perturbations in models with decaying fields. 

The main simplifying approximations we made were the neglection of perturbations in radiation 
and the use of a sudden decay approximation, requiring the inclusion of a window function ()65p 
in the matching conditions to induce a suppression of the ringing patterns on small scales. Any 
surviving ringing patterns on small scales should be taken with caution, since we ignored terms 
proportional to (fc/'H)^ in deriving the matching conditions for perturbations.^^ 

We focused on a simple multi-field model, that, although motivated by moduli dynamics on the 
landscape, should primarily be seen as a toy model to test the tools developed so far to deal with 

Note that in [sj the terms of order 0{{k/T-Cf), which we keep, are ignored. 
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staggered inflation. 



V. CONCLUSION 

We computed numerically the scalar power-spectrum in a simple multi-field inflationary setup 
that contains decaying fields during inflation, varying the number of fields from a few to thousands, 
motivated by multi-field models in string theory. For low field numbers, we find distinguishab^ 



33|, 



patterns in the power-spectrum caused by individual decays, in line with the analytic results of 
where a single decay has been investigated in a concrete string-derived model. Once many fields 
decay durin g an y g iven Hubble time, we recover the results of the analytic treatment of staggered 



inflation in [38|, |40(]: the superposition of many patterns yields a smooth spectrum that is well 
described by a power-law; the scalar spectral index is not only set by the slow roll parameters, but 
also by the ratio of the decay rate to the Hubble parameter which, in the model under consideration, 
is the dominant contribution. 

The numerical confirmation of the analytic results in the large/small M limit, as well as the 
development of a reliable numerical code to deal with any number of fields are the main results of 
this paper. 

A shortcoming of our current approach is the neglection of perturbations in the radiation bath 
that is produced whenever a field decays. This made an overall renormalization of the Bogoliubov 
coefficients a necessity whenever a field decayed and leads to order one deviations for the scalar 
spectral index as compared to the analytic results of staggered infiation. We do not expect qualita- 
tive differences if these perturbations were kept. We also ignored terms of order k^/T-l^ and worked 
in a sudden decay approximation, which prevents a suppression of ringing patterns on scales deep 
within the horizon. To ameliorate this unphysical effect, we introduced a window-function into 
the matching conditions by hand, rendering modes deep inside the horizon oblivious to the fields' 
decay. 

We plan to use the code developed for this paper in an investigation of the extended KKLMMT 
setup in [33], where we plan to relax these approximations. 
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Appendix A: Code Implementation 

To evolve the scale factor computationally the Friedmann equation and Klein-Gordon field 
equation were rewritten as coupled first-order differential equations. These modified equations 
take the form 

dVr 

= -7J-^-3e^/, (Al) 



^ = ^^(EQ^K^.)+pJ- (A2) 

We track the natural logarithm of the scale factor, ^ = In a, which is the number of e-folds since 
the beginning of the simulation. 

Similarly, the second-order differential equation for the Sasaki-Mukhanov variable (j43p is cou- 
pled to ()A2p and can be decomposed into two coupled first-order differential equations. The coupled 
first-order Sasaki-Mukhanov equations are 



N 



(A3) 



J 

Since the Sasaki-Mukhanov variables are complex, we track both the real and imaginary compo- 
nents. At the beginning of the simulation we chose a set of ki . . .ki where i w 100 — 500. This set 
of discrete Qi{ki) are initialized in the Bunch-Davies vacuum, (|40p and the Qi{ki) are identical for 
each /. All modes are tracked analytically until the physical wavelength a/k is some fraction /3 of 
the Hubble length, i.e. until k = aH/f3; for the simulations in this paper we take f3 = exp(— 4). At 
this point we evolve the Qj numerically according to (jASp . 

The equations are evolved using the fifth-order Cash-Karp Runge-Kutta method of integration 



64l |. The primary advantage of this method is its adaptive step size. At each time step the program 
computes values via fifth-order integration and an embedded fourth-order integration. The fastest 
oscillating mode in play determines the adjustable time-step. All other modes in play oscillate at 
a lower frequency and are therefore not as sensitive to the timestep. 
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